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Summary 

This computational study is an exploration of 
nonreacting high-speed mixing layers. The compu- 
tations are made with a code based on an accurate 
higher order algorithm and with sufficient grid points 
to resolve all relevant scales. In each case, the free- 
stream disturbance is introduced, and the calcula- 
tions are run until a statistical steady state is reached 
(2 or 3 sweeps, based on axial distance and convec- 
tive speed). Then, 400 to 600 time samples of the 
flow, at equal time intervals, are generated to obtain 
statistical properties of the flow. 

The studies are conducted at two convective Mach 
numbers, three free-stream turbulence intensities, 
three Reynolds numbers, and two types of initial 
profiles — hyperbolic tangent (tanh) and boundary 
layer. These profiles are used to interpret the effects 
of corresponding parameters. 

The boundary- layer profile leads to predictions of 
transition processes much better than the tanh pro- 
file at disturbance intensities comparable to those in 
wind tunnels. The transition Reynolds number pre- 
dicted is about 0.18 x 10 6 and compares well with 
experimental data. The spacing between the vorti- 
cal structure normalized by the boundary-layer thick- 
ness is about 3.5 (boundary-layer case) and 10 (tanh 
case), whereas the observed values range from 1.5 
to 2.5. The deduction of the growth-rate data from 
the plots of growth versus the streamwise coordinate 
is, in part, affected by the nature of free-stream dis- 
turbance introduced and the regions of slight nega- 
tive growth. Nevertheless, careful extraction of the 
growth-rate information and comparison with the 
experimental data of Papamoschou show that this 
comparison is moderately good. The study of the 
location of the disturbance shows that the growth 
of the layer is suppressed if the shear-layer region 
is excluded from the disturbance field. Calculations 
with the boundary- layer case with no imposed distur- 
bance show that the shear layer selectively amplifies 
the noise at a Strouhal number (based on momen- 
tum thickness) of 0.007. This result can be used in 
practice to excite the shear layer toward early tran- 
sition. Studies of dependence on Reynolds number 
show that the effects are negligible, and this may be 
caused by the dominant inviscid structure of the flow 
field. 

Introduction 

High-speed mixing layers have received significant 
attention in the past 5 years because of the poten- 
tial applications in hypersonic plane scramjet com- 
bustors, where the combustion process is supersonic. 


Recent experimental results (ref. 1) indicate that the 
growth rates of high-speed layers in relation to those 
of incompressible shear layers are much smaller— by 
a factor of 4 to 5. This reduction implies that much 
longer lengths are needed to achieve the same level of 
mixing and combustion efficiency. This longer length 
is obviously undesirable from a technological point of 
view. Consequently, considerable attention is being 
paid to the mixing phenomenon in supersonic mixing 
layers — both nonreacting and chemically reacting. 

Studies of stability of mixing layers (refs. 2 
and 3) and direct simulation of temporally and spa- 
tially developing mixing layers (refs. 4 and 5) are 
the theoretical and computational approaches made 
to understanding the behavior of the mixing layers. 
Experimental studies of the high-speed mixing layers 
are difficult, and there are not many such studies. 
These mixing layers have not been sufficiently ex- 
plored for the distribution of mean and fluctuating 
quantities with nonintrusive diagnostic tools. The 
only reliable sources of information available are the 
growth rates measured from pitot traverses, Schlieren 
pictures, and mean velocity measurements. Some ef- 
forts at modeling the mixing process have been made, 
but a more precise understanding of the flow must be 
obtained before the modeling efforts are considered 
valid. 

In light of this scenario, direct numerical simula- 
tion of the shear layers becomes an important tool 
for understanding the mixing process. The classical 
question of resolution of scales at high Reynolds num- 
bers limits the exploration to relatively low Reynolds 
numbers («1 xlO 4 ) based on shear-layer thickness 
and mean velocity at the inflow plane). But this 
range is not much below the range of practical rel- 
evance or importance. Thus, it should be possible 
to perform these calculations and to make deduc- 
tions of relevance and some importance to the be- 
havior of mixing layers. Initially, studies will be 
limited to two-dimensional simulations. While there 
is the capability to make three-dimensional simula- 
tions, it is important to understand two-dimensional 
simulations before embarking on three-dimensional 
calculations. In recent times, a few groups (refs. 4 
and 6) have made direct numerical simulations of 
shear layers, and much has been uncovered. The 
present contribution is complementary in some as- 
pects and presents new results on many aspects. 

Symbols 

a speed of sound 

b.l. boundary layer 

c constant 



fi mass fraction of species i 

h enthalpy 

k turbulent kinetic energy 

M Mach number 

M r convective Mach number 

p pressure 

R Reynolds number 


Previous Studies 

Experimental 

Table 1 contains a list of references and some 
brief details of the experimental work discussed in 
them. This table has been constructed to determine 
which experiments can be used for comparison with 
the results of the simulation. Figure 1 shows the 
experimental design used by these authors to obtain 
their results. 
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transition Reynolds number 
root mean square 
Strouhal number 
temperature 
time 

hyperbolic tangent 

primary stream velocity 

secondary stream velocity 

stream wise velocity 

con vec t i ve ve loci ty 

disturbance in streamwise direction 

streamwise mean velocity 

transverse velocity 

disturbance in cross-stream direction 

streamwise coordinate 

maximum extent of x field 

transition distance 

transverse coordinate 

maximum extent of y field 

ratio of specific heats 

mixing- layer thickness 

growth rate 

Kolmogorov scale 

momentum defect thickness 

laminar viscosity 

density 

vorticity 


Subscripts: 
p primary 

5 secondary 


Reference 7 is a critical summary of the data of 
many early studies on both subsonic and supersonic 
shear layers. The data of references 8 and 9 and their 
nondimensional growth rate as a function of Mach 
number are important. The growth rate is normal- 
ized by that under incompressible conditions. It is 
apparent from reference 10 that the incompressible 
mixing-layer growth rate is uncertain to a significant 
extent (a factor of 1.8). Therefore, the quality of 
the results of nondimensional growth rate versus M 
is not assured, even though the primary result that 
mixing-layer growth rate decreases with convective 
Mach number M r is unquestioned. 

Experiments have been conducted on mixing lay- 
ers (ref. 11) by using a splitter-plate configuration 
(fig. 1(a)). Measurements of mean velocity profile 
arid transition have been made for two configurations 
at M c = 0.325 and 0.515 (table 1), and the pressures 
across the mixing layer are not the same. In one 
case, they are different by a factor of 1.3. Transi- 
tion is measured in the following ways: (1) point of 
change in slope of the plot of growth versus distance, 
(2) examination of Schlieren pictures, and (3) point 
of decay of the spectral intensity of the dominant 
frequencies in the free shear layer as measured by 
a hot-film probe. The transition locations, as mea- 
sured by these criteria, are not necessarily the same-- 
they vary by about 35 to 45 mm. The transition 
distance normalized by initial momentum thickness 
is about 165 to 300#, depending on the flow condi- 
tions and the transition Reynolds number based on 
mean properties (and mean velocity), which is about 
2.0 x 10 5 ± 0.3 x 10 5 . 

The experiments of reference 12 are instructive in 
many ways. Though the study is limited to a sin- 
gle speed and the geometry used is not pertinent to 
the two-stream mixing considered here, it gives an 
appreciation of the free-stream turbulence and the 
peak turbulence intensity inside the shear layer. The 
measured free-stream turbulence is about 0.2 percent 
for velocity fluctuations and 0.1 percent for pressure 
fluctuations. The peak turbulence intensity increases 
to 6 percent in the downstream region of the shear 
layer. This study also presents the velocity spectra 
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that can be used for detailed comparison. The appa- 
ratus used by these spectra is shown in figure 1(c). 
The injection from the bottom wall helps balance the 
pressure across the shear layer. 

Papamoschou’s work (refs. 1 and 13 to 15) is by 
far the most referenced work, essentially because of 
a systematic effort to evolve appropriate parameters 
for compressibility effects and to determine the ef- 
fect of compressibility on mixing (growth rate of the 
mixing layer). There were previous studies (both ex- 
perimental and linear stability) on the subject, but 
the status of understanding was not consolidated un- 
til his work clarified issues and reinforced the impli- 
cations of previous work on stability (ref. 2). The 
principal result he obtained was the variation of the 
nondimensional growth rate with convective Mach 
number. Following the work of Roshko and others, 
the dynamics of the mixing layer is argued to be af- 
fected by the speed of large-scale structures relative 
to the speed of one of the streams. The ratio of this 
speed to the average acoustic speed is the convec- 
tive Mach number M c . The growth rate normalized 
by that for incompressible conditions decreases by 
a factor of about 4 to 5 when M c increases to 2. 
In these experiments, Papamoschou used nitrogen, 
argon, and helium as test gases in a splitter-plate 
configuration (fig. 1(b)). Understandably, the test 
durations were short, typically 1 to 2 sec, and the 
measurements consisted of pitot surveys across the 
mixing layer. The free-stream turbulence effects are 
only indirectly inferred from wall pressure measure- 
ments. One of the other details measured, but not 
published, is the difference in the static pressure of 
the two streams. This difference may well be small, 
but it adds to the wave processes and, in a disturbed 
environment, may affect transition. 

Some preliminary measurements of a shear layer 
of a single fluid with two streams across a splitter 
plate are reported in reference 16. A single measure- 
ment of growth rate and velocity profile has been 
reported, as well as a mean turbulence velocity of 
5 percent in the mixing layer. Reference 16 also in- 
dicates from the measurement of laser doppler ve- 
locimetry (LDV) that the wave reflections from the 
walls can cause fluctuations of flow velocities as high 
as 4 percent. 

Some experiments of supersonic shear layer have 
been presented in reference 17. In this setup, the 
shear layer begins from a backward-facing step, and 
the flow structure in the downstream zone is altered 
by inserting a ramp (fig. 1(d)). The ramp provides 
a recirculation zone at near-zero speeds and helps to 
maintain a smooth shear layer. The growth rates 
from these experiments, with the incoming stream 


deliberately made turbulent, are plotted along with 
the previous data to reinstate the conclusion regard- 
ing the reduced growth rate at supersonic conditions 
(in terms of M c ). Transition experiments have been 
conducted recently on a similar setup (ref. 18), but 
the test section was made disturbance free to a large 
extent. Measurements of transition Reynolds num- 
ber show values from 3.63 x 10 5 to 5.3 x 10°. 

It is clear from the preceding discussion that 
(1) the various experimenters have used different ap- 
paratus for creating a mixing layer, (2) the measure- 
ments of the quality of the incoming stream are not 
complete in most cases, and (3) the mixing layers are 
subject to disturbance, acoustic and otherwise, at a 
magnitude upward of 0.2 percent, and a pressure dif- 
ferential across the layer of 1 to 1.3 in some instances. 
Yet, the growth rates of the layer from the work of 
most investigators, possibly in the turbulent range, 
are put on the same plot, and evidence of these be- 
ing nearly the same qualitatively and quantitatively 
is presented (ref. 13). If these are true, then the 
asymptotic growth rate is relatively independent of 
tunnel characteristics and the disturbance environ- 
ment. This aspect needs confirmation through sys- 
tematic experimental study in which the incoming 
disturbance field is measured along with the features 
related to growth rate. 

Computational 

The efforts to model the flow, particularly to ex- 
plain the reduced growth rate, are recent. These 
efforts have been related to stability, direct sim- 
ulation, and a correlation of both. References 3 
and 5 show that, at high speeds, the peak ampli- 
tudes of the unstable waves decrease with increas- 
ing M c , Also, at high values of M c , the three- 
dimensional oblique disturbances tend to grow more 
than the two-dimensional waves. This reduced am- 
plification of the disturbances may be responsible for 
the delayed growth (in terms of distance from the 
edge of the splitter plate) and, in general, the inabil- 
ity of the flow to sustain the growth of disturbances. 

Direct simulation of mixing layers has been con- 
ducted in reference 6 for supersonic mixing layers 
and in reference 19 for incompressible mixing lay- 
ers. A high-order finite- difference scheme is coupled 
with fine-resolution grids to capture the important 
scales of the transport within the mixing layer in ref- 
erence 6. The computations are aimed at establishing 
the importance of M c and elucidating the causes for 
reduced mixing at high Mach numbers. The com- 
putations use tanh profiles and sufficiently large dis- 
turbances on the cross-stream velocity at the inflow 
station. The inflow disturbance magnitude is as high 
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as 5 percent. While some conclusions are not af- 
fected by this unnaturally large forcing, others are. 
For instance, the mixing-layer thickness versus axial 
distance shows a linear growth that tends, asymptot- 
ically, to zero growth. This result has no parallel in 
reality and is left unexplained. It is concluded from 
the present simulation, as well as from reference 19, 
that this result could be one of the effects of forcing. 

The results for large Mach numbers indicate the 
presence of thin shocks called eddy shocklets. The re- 
duced mixing at high Mach numbers has been related 
to the change in balance of terms involving baroclinic 
and compressibility terms in the velocity equation. 
Arguments of this nature were made in the temporal 
simulations of reference 20; these arguments indicate 
that baroclinic torque resulting from compressibility 
causes distribution of vorticity; the diffusion of the 
vorticity over a wider region results in weaker vortex 
roll-up and slower growth. 

Reference 19 contains numerical simulations on 
incompressible mixing layers. The spread rate has 
nonmonotonic behavior with distance for specific 
cases, and this is related to the phase relationships 
between forcing and the pairing processes. This 
behavior is also related to the change in sign of 
the Reynolds number stress. Specifically, in regions 
where the spread rate decreases or remains nearly 
constant with distance, the energy is transferred from 
large-scale structures to the mean flow. 

Reference 10, a review of free shear-layer mixing, 
addresses several relevant issues. Many aspects of 
this paper are for subsonic flows. The incompressible 
shear-layer growth rate is expressed as 


where c$ is a constant, possibly independent of veloc- 
ity ratio r and/or the density ratio s . Based on cu- 
mulative experimental evidence, it is concluded that 

varies between 0.25 and 0.45. The large difference 
is still unaccounted for. Experimental evidence from 
incompressible flows is presented to indicate that the 
growth rate is a strong function of Reynolds number 
in a mixing layer. 

The growth-rate dependence on Reynolds num- 
ber seems to be a subject of some controversy. 
Papamoschou discusses this question in reference 13; 
he does not take a firm position, but indicates that 
he ignores the effect of Reynolds number. The reason 
that the large dependence on Reynolds number ap- 
pears difficult to accept is that the mixing-layer phe- 
nomenon, particularly in the early part, is largely 


inviscid. Hence, one would expect weak Reynolds 
number dependence at best. 

The present work addresses the question of 
Reynolds number dependence by changing the 
Reynolds number through alterations of pressure 
rather than velocity. This approach has been dis- 
cussed by Dimotakis (ref. 10), and it is one of the 
cleanest ways of changing the Reynolds number with- 
out altering other aspects of the flow. For example, 
changes of velocity also alter the initial boundary- 
layer profiles and have influences which have not been 
evaluated until now. One of the aims of this work is 
to treat this aspect as well. 

Outline of Present Work 

The present work is a direct simulation of mixing 
layers with different reactive fluids in the two adjoin- 
ing streams. They are chosen so as to form a set for 
reactive studies as well. The first part of the work 
discusses the methodology for analyzing and inter- 
preting the results of computation. The questions 
concerning grid resolution and boundary effects are 
discussed subsequently. The results of effects of free- 
stream disturbance on the growth of the shear layer 
at two different values of M c are then brought out; 
both tanh and boundary layer, like initial profiles, 
are discussed. The various aspects — growth rates of 
vorticity, velocity, mass fraction, and density; vor- 
ticity plots and variation of turbulence quantities 
through the flow field; and time and space spectra 
of the fluctuations are examined to obtain insight 
into the behavior of the mixing layer. The cases con- 
sidered in the present work are shown in table 2. 
The choice of the two cases is aimed at making cal- 
culations, one at low M c and another at a relatively 
high value. For the two cases considered, M c — 0.38 
and 0.76. A fuel and oxidizer system was chosen in- 
stead of single fluid system to enable computations 
with reaction as well (not reported herein). The tem- 
perature of both streams is held the same, so as not 
to additionally introduce a parameter through tem- 
perature. The temperature for the present test is 
2000 K, as it is representative of inflow conditions 
in a supersonic combustor of a hypersonic cruise ve- 
hicle. Table 2 shows that the specific-heat ratio 7 
is nearly the same for both streams, = 1.0133. 
This similarity is because of the high temperature of 
the stream and because the molecules involved are 
diatomic in nature. 

The inflow-plane pressure is the same on both 
sides of the mixing layer at 0.101325 MPa (1 atm). 
Changes in Reynolds number are made by changing 
the pressure above and below the nominal value by a 
factor of 2. The convective properties are calculated 
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from 


and 


( m p - m *v 15) 
te)*( i+ v^) 


Uc = 


yJ~PpUp ~h yj PsUs 

{y/Pp + y/Pi) 


(2) 

(3) 


These are the definitions presented in reference 14. 
They use the geometric average of the values of M c 
on the primary p and secondary s sides. The primary 
side is defined as the one with the higher Mach 
number, no matter what the speed. In equation (2), 
7 has little impact. If we ignore the 7 terms, then 
equation (2) reduces to 


= y/pIMp (4) 

\fPp + sfFs 

As discussed previously, Papamoschou s definition of 
M c (ref. 13) is 


M c 


Upj^CJs 

dp + d s 


(5) 


Equation (2) reduces to the above result if we recog- 
nize that a = sjip/p, and p and 7 are the same for 
the two streams. The initial profiles are chosen as a 
tanh profile as follows: 


u m = 1 [(Uoc + U- 00) + (E/00 - U- 00) tanh ky] (6) 

The constant k is chosen as 1800 m \ so that 
the effective boundary-layer thickness, based on a 
99-percent free-stream velocity criterion, is about 
2 mm. For other numerical experiments, boundary- 
layer profiles are given by 

u m = t/oc[l -exp(— fciy)] (y>ym/ 2) (7) 

and 

u m = U -00 [1 - ex P {-hy)\ ( y < W 2 ) ( 8 ) 

The constant Jt, is chosen as 4000 m^ 1 . The effective 
boundary-layer thickness on each side is about 1 mm. 
This profile is intended to simulate the normally 
obtained experimental conditions. 

Table 3 shows the nondimensional parameters 
relevant to the cases considered. Several instructive 
features can be derived from tables 2 and 3. The 
velocity of the airstream is lower than that of the fuel 
for M c = 0.38, but is much higher for M c = 0.76. The 


momentum ratios shown in table 3 indicate that the 
momenta are balanced for M c = 0.38 but are greatly 
in favor of the airstream for M c = 0.76. The shear 
layer therefore remains roughly in alignment with the 
central axis in the first case, but bends over toward 
the fuel side in the second case. The implications 
for the computations are that the grid has to be 
well resolved, even in the outer regions, if the details 
are to be captured in the second case, or a twisted 
stretched grid must be used. In fact, this problem 
can be so serious that reactive flow computations do 
give considerable difficulty in the second case unless 
one of these strategies is adopted. 

The second part of the work concerns the layer 
thicknesses. The thicknesses, based on 99-percent 
free-stream velocity, are about 1 mm on each side 
of the splitter plate. The momentum defect thick- 
nesses 6 for the tanh and boundary- layer profiles are 
very different from each other. In the tanh case, it is 
about 1.54 mm, and in the boundary-layer case, it is 
about 0.102 mm. This large difference implies that 
for the same streamwise grid length (e.g., 100 mm), 
only about 64 momentum thicknesses can be cov- 
ered with the tanh profile, but about 1000 momen- 
tum thicknesses can be covered with the boundary- 
layer profile. For small free-stream disturbance levels 
consistent with experiments, no transition should be 
expected for the tanh case, but significant transition 
should be expected for the boundary-layer case. This 
situation can be altered if a free-stream condition of 
high turbulence intensity is imposed. If the transi- 
tion is caused early enough, the asymptotic growth 
rate from this calculation would be representative 
of turbulent boundary-layer growth. The Reynolds 
number based on boundary- layer thickness is about 
3.5 x 10 3 to 1.2 x 10 4 for the two Mach number cases. 
Based on the convective speed, momentum thickness, 
average density, and viscosity, then R = 360 to 550. 
If the reference speed is changed to the difference 
in speeds between the two streams, then R = 250 
to 400. This last quantity is brought out because the 
difference in velocities is used by some workers for 
normalization purposes (ref. 6). These values are low 
enough that the direct simulation approach with the 
current-day computational aids is expected to lead 
to realistic results. 

For disturbances to be introduced with the inflow 
profiles, several numerical experiments have been 
conducted. The disturbance function is chosen as 

3 

u d = (Uoo - ( 9 ) 

2—1 
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Vd = (Uoo - U- op) dj sin Ujt (10) 

i=l 

The velocity variation with y at x = 0 is then given 
by u = 4- lid and i; = In most cases, 

f-2 = C 3 — 0, = 100 kHz, ci = 0.04, d\ — 0.04, 

d 2 = 0.02, d 3 = 0.01, a ?2 = 200 kHz, and 

= 300 kHz. 

In a few preliminary runs made with no inflow dis- 
turbance, the layer grew very little to about 100 mm. 
There was only a trace of growth toward the end. 
The time-sampled data at this plane were subjected 
to spectral analysis, and the frequencies, with rel- 
atively large amplitudes, were chosen and used in 
equations (9) and (10). The choice of the amplitudes 
is arbitrary and is made to ensure that the rms val- 
ues of the fluctuations are of a required magnitude. 
The above disturbances (u ( p v^) are introduced ei- 
ther into a specific region of about 4 mm or over 
the entire region. In some experiments, the distur- 
bances are restricted to the fuel or oxidizer region 
only. These calculations were performed to observe 
the sensitivity of the mixing layer to the location of 
the disturbances. 

Code and Algorithms 

The code used in the present calculations is the 
SPARK combustion code developed at Langley Re- 
search Center over the past 4 years. It solves 
the Navier-Stokes equations, including energy and 
species conservation. The models chosen to eval- 
uate viscosity, thermal conduction, and mass diffu- 
sion allow for a mixture of perfect gases. In the first 
version, the code used a second-order spatially and 
temporally accurate, two-step explicit MacCormack 
scheme. (See ref. 21.) Subsequently, it has been mod- 
ified to include a variety of higher order spatial al- 
gorithms, including linear and compact algorithms 
(fourth and sixth order) and various upwind algo- 
rithms (third and fifth order). A brief description of 
the algorithms is provided in the appendix. Based 
on the supersonic streamwise characteristics of the 
mixing layer, a third-order, upwind-biased algorithm 
is used for the streamwise direction. A fourth-order, 
central-difference algorithm (ref. 22) is used in the 
cross-stream direction. The temporal accuracy is still 
second order. This choice represents a compromise 
between the accuracy of higher order numerical al- 
gorithms and the robustness and efficiency of lower 
order methods. 

Boundary Conditions 

The problem is cast in the x-y coordinates, such 
that the flow is oriented along the x-axis. The two 


streams extend from y — 0 to y rn / 2 and y m / 2 to 
2/m, where y m is the maximum extent of y. In the 
streamwise directions, the boundary conditions at 
x = 0 are the set values of u, v , p, T, and f u 
since the inflow is supersonic. The fluctuations, when 
introduced, are only in u and v. In a few experiments, 
fluctuations in p are also set. At the supersonic 
outflow plane x = x m , first-order extrapolation of 
the primitive variables is used. 

At y = 0 and y m , the gradient ^ r °P^ Tt y) _ q j s 
set. To ensure that the boundaries are far enough 
for this condition to be valid, the outer regions 
are chosen as 15 times the initial boundary- layer 
thickness. Therefore, y m = 30 rum. A few numerical 
studies have been conducted with y m = 50 mm. In 
one case at M c = 0.76, doubts were raised as to 
whether the results were affected by the boundary. 
In this case, the boundary line is set at y m = 100 mm. 
In all cases, the results obtained at y rn — 30 rnm were 
reproducible arid, to this extent, seem satisfactory. 

Choice of Grids 

It is intended that most scales of importance be 
captured in the calculations. Numerical resolution of 
the large-scale vortical structures presents no prob- 
lem. However, the relatively fine scales at which dis- 
sipation through viscosity occurs need fine grids for 
capture. The Kolmogorov scale, which describes the 
fine scales, is related to the large scale (integral scale) 
as follows: 

7/ 3 

j=CR’* (11) 

where R is the Reynolds number based on a charac- 
teristic thickness. The Reynolds number of the flow 
with the initial boundary- layer thickness as the rele- 
vant scale varies from 3.5 x 10 3 to 1.2 x 10 4 for the 
problems treated. In the downstream region, where 
the shear-layer thickness is larger, Reynolds number 
goes up to 1.5 x 10 4 . The constant C is to be es- 
timated. The balance of dissipation and production 
is used in reference 23 to estimate the constant as 
typically 2.5 for incompressible flows. Recognizing 
the fact that in compressible flows the disturbances 
tend to have much lower growth rates — typically 2 
to 3 times smaller than in incompressible flows — and 
that the asymptotic structure retains this effect, the 
constant is between 5 and 7. For a typical extent 
of growth estimated from known data — 6 =10 mm 
and R = 1.0 x 10 4 — 77 = 0.056 mm. The smallest 
size of the computational domain is 0.025 mm. The 
stretching of the grid makes the grid size increase 
to 0.4 mm at a distance of 10 mm from the central 
line. This increase does not ensure that all the scales 
will be captured, as calculated in reference 17, but it 
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does ensure that most scales of importance are cap- 
tured. Since this is an important issue, the justifica- 
tion for the choice of the step size is not limited to 
this argument. The spectra obtained from two differ- 
ent fine grids are compared to provide the necessary 
justification. 

The number of grid points used in the y direction 
varies from 101 to 176. Many results obtained here 
are for 151 grid points over a region of 30 mm. A 
grid compression factor of 10.0 allows the region of 
10 mm to be covered with about 95 grid points. 
The x direction is typically embedded with 0.5-mm 
equispaced grids. Some numerical experiments with 
a grid-point spacing of 0.25 mm were also obtained 
to ensure that the results are grid-independent for 
a large number of fluctuating parameters that are 
examined. The total number of grid points is 201, 
251, or 301, depending on the region examined 50, 
100, or 125 mm, respectively. In one case, a 200-mm 
region was treated, because the growth of vortical 
structures was virtually absent up to 125 mm. The 
minimum step size chosen should, again, be able to 
capture the scales. The vortical scales in the axial 
direction are 3 6 to 46. These scales are covered by 
60 to 80 points, which is good in most cases where 
the flow is in the transitioning stage, but marginal in 
the boundary-layer case in a zone where turbulence 
has set in significantly. 

Methodology and Tools of Examination 

It is important to understand the flow behavior in 
each of the cases computed. Since there is a variation 
of the flow behavior in both space and time, there 
is a need to analyze both. The spatial variation 
is understood largely by contour plots of vorticity, 
pressure, and fuel mass fraction. The presence of 
shocks is best understood in terms of a quantity 
called the “shock function based on pressure.’ 1 Other 
plots, including the v velocity- the divergence of 
velocity- are usually less revealing with regard to 
the presence of sharp changes of pressure in the 
flow. However, this presentation is restricted to the 
contour plots of pressure. 

The speed of vortical structures (expected to be 
near the convective speed) is obtained from pressure 
versus x plots at some value of y near the center. 
These plots are made at a few times that are suffi- 
ciently spaced apart. The rate of movement of the 
point of peak pressure along the x-axis is the convec- 
tive speed. 

To determine a statistical time behavior, it is nec- 
essary to time sample the flow, preferably at equal 
time intervals. The SPARK2D code calculates the 


time step based on a numerical stability criterion 
(CFL). It varies significantly in the early part of the 
calculation but settles down to a value with a fluc- 
tuation of about 5 percent. To take these features 
into account and the fact that the flow must attain a 
statistical steady state before sampling is performed, 
the code is run for each case for a duration of about 
three convective sweeps of the flow. Each sweep takes 
a time given by x m /U ( \ this is about 50 msec for the 
M r = 0.38 case and 30 msec for the M r = 0.76 case. 
The time step is typically 0.005 //sec; therefore, it 
takes 20000 to 30 000 time steps before statistical 
steady state is achieved. A total of 400 to 600 time 
samples of the flow field involving p, u, i\ T, p, //, 
and vorticity, in the most important y region and 
at specific x stations, are stored at equal time inter- 
vals. These samples are subsequently analyzed by 
a separate statistical package specifically developed 
for this purpose. The results from this package? in- 
clude several quantities: the mean and root mean 
square of fluctuations of all the parameters, Reynolds 
stress and other correlations of velocity with tem- 
perature, mass fractions and pressure, average thick- 
ness estimated by various means, probability density 
function of the passive scaler (mass fraction of H 2 ) 
and other parameters as needed, time spectra, and 
x spectra (wave-number plot) for the variables de- 
sired. The intention was to obtain the spectra in the 
y direction, but this objective was impeded by the 
problem of calculating the error-free integrals of the 
variables known at specific points (unequal intervals 
in y) with harmonic functions. These calculations 
produced a significant amount of high-frequency 
noise, even though a simple quadrature? rule was used. 
The x spectra obtained from data of equal intervals 
were sufficient to draw inferences. The shear-layer 
thickness was obtained for u , vorticity, density, and 
H *2 mass fraction. Of these, only vorticity tends to 
zero at y — ± 00 , and the others tend to nonzero 
finite values. Therefore, the thicknesses are defined 
by 

c {V^i-U-oc) /Uo\ 
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with similar definitions for density and mass fraction 
of H 2 , with the velocity u replaced by the correspond- 
ing quantities. For vorticity thickness, the definition 
is given by 
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(13) 


where f) is vorticity. 

Defining 6 as in equations (12) and (13) turned 
out to be reasonable for several cases. However, in 
the case of the boundary layer and in some cases of 
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tanh profile, the preceding definitions gave what ap- 
peared to be inconsistent results. In the transitional 
regime (a regime where the vortical structures begin 
to interact and lead to other structures of smaller and 
larger sizes), the development of the profiles was such 
that the peak gradients of mean quantities did not 
possess a monotonically decreasing trend. A small 
kink in the profile, as can be seen in some of the plots 
of mean quantities described subsequently, upset the 
inonotonicity of the gradients. As such, a spurious 
decrease of the growth of the layer was indicated, 
even though an examination of the plot of the mean 
profile did not show such a behavior. For this reason, 
a thickness based on a 98.5-percent free-stream value 
was obtained and presented. 

Results and Discussion 

Hyperbolic- Tangent Profile 

Figure 2 shows a typical result of a structure 
superimposed on the grid and the grid distribution. 
The maximum size of the vertical structure is about 
10 mm, the minimum grid size is 0.025 mm, and the 
grid size at the outer edge of the structure is 0.4 mm. 

Figure 3 contains the plots of growth rates with 
axial distance for the nominal case ( M c — 0.38) along 
with a finer y grid and a finer x grid. All the cases 
are for a region 100 mm (in x) by 30 mm (in y). The 
nominal case is for 201 by 101 grid points. The finer 
y- grid case has 151 grid points. The finer x-grid case 
has 301 grid points. The growth rates are obtained 
by averaging 250 time samples. The differences in 
the results of figure 3 are the results of averaging and 
predictions of grid resolution at the levels employed 
for the nominal case. 

Figure 4 contains similar plots to examine the 
influence of the boundary. The y region of 30 mm 
is extended to 50 mm. Again, the differences seen 
are simply the result of statistical averaging, and 
the boundary of y = 30 mm should be satisfactory. 
The question of the outer y boundary being kept at 
Vm = 30 mm (15 mm away from the center of the 
shear layer) becomes more serious for the M c = 0.76 
case, because, in this case, the changes in pressure 
are much larger than in the M c = 0.38 case, and 
boundary effects on the radiation of these pressure 
waves could be more significant. 

Figure 5 contains the contour plots of vorticity 
and pressure for the small domain (y m = 30 mm 
with the mixing layer centered at 15 mm) and the ex- 
tended domain (y m = 100 mm with the mixing layer 
centered at 50 mm). These plots use the same scale 
and same range of vorticity at about the same time, 
though not identical. The fact that the vorticity plots 


look identical may be misleading. Also, the fact that 
the pressure plots retain nearly identical behavior is 
a better indication that the outer-boundary set at 
Vm = 30 mm is reasonable. Other data involving 
means and fluctuations in the central 10 mm region 
show differences less than 3 percent not very much 
larger than errors of statistical averaging, which are 
about 2 percent. The verification of grid resolution 
and boundary effect based on growth rates, means 
of flow variables, and contour plots of pressure was 
thought to be more complete with the examination 
of time and wave spectra. 

Figures 6 and 7 refer to time and space spectra 
for u and v shown at the initial and final stations. 
The plot at the initial station shows the frequencies 
used to excite the layer in u and v. The spectrum 
approaches a function with a sharper peak by us- 
ing a larger number of samples. The differences in 
spectra seen in other cases are the usual statistical 
differences. The spectra show more than a two-order 
range of frequencies from about 12 kHz to 1.5 MHz. 
The peak energy containing eddies have a frequency 
between 40 and 100 kHz. The quantity shown as 
amp is the root of kinetic energy of u fluctuations 
in the specific frequency range. The dominance of 
the initial disturbance at a frequency of 100 kHz is 
seen 50 mm downstream, but generation of both sub- 
harmonic and higher harmonic frequencies is taking 
place. At 100 mm downstream (600), the amplitude 
of the initial disturbance has come down from 100 
to about 50 m/sec, and the other lower frequencies 
have increased in amplitude substantially (i.e., from 
about 5 m/sec to as much as 20 to 50 m/sec). Also, 
the higher frequency spectrum is populated. 

Figure 8 contains the velocity-vector plots, con- 
tours of vorticity, and pressure and mass fraction of 
H 2 for the tanh profile case at M c = 0.38. These 
results are included to illustrate the relationship be- 
tween various quantities. The mass fraction of hy- 
drogen follows the vorticity plot, which shows that 
mixing is related to the vortical structures. The pres- 
sure plot shows that the pressure is lowest near the 
center of the vortical structures and peaks between 
the structures. 

The influence of intensity of fluctuations on the 
growth of vortical structures is shown in figure 9. 
With no inflow disturbances, the growth is negligible. 
The viscous entrainment in the laminar range is, 
indeed, extremely low. With an increased level of 
disturbances at 0.92, 1.84, and 3.68 percent (rms 
value), the instability begins at smaller distances, 
and the growth of vortical structures is stronger. 
If the growth of the structures was associated with 
the transition process, then it could be inferred that 
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the process has occurred at about 80 to 90 mm at 
0.9 percent rms disturbance. For realistic transition 
distances of the order of 40 to 50 mm, the intensity of 
fluctuation needed is of the order of 3 percent. This 
is very large compared with the kind of free- stream 
disturbances expected in wind tunnels or free flight. 
This feature is strongly related to the nature of the 
initial profile. 

Figures 10 and 11 contain the corresponding plots 
for M c = 0.76. These plots are included for compari- 
son with those for M c = 0.38. The growth rate of the 
vortical structures for 3.6 percent fluctuation seems 
as large as for the M r = 0.38 case. 

Figure 12 shows the effect of the location of the 
disturbance on the growth of the shear layer. Dis- 
turbances in the free stream away from the shear 
layer seem to have little influence on the growth of 
the layer. Inclusion of the shear layer into the dis- 
turbance field makes all the difference, as seen from 
the third and fourth plots. Though a careful exam- 
ination of these plots shows that the disturbance in 
the air side (lower) leads to slightly larger structures, 
a generalization is not possible with the limited re- 
sults available. Disturbances set around a region 
that does not include the shear layer are not suffi- 
cient to realize the effects of the disturbance field. 
It is necessary to include the shear layer and a re- 
gion of about four times the layer thickness to obtain 
the full effect of the disturbances. This result may 
have important practical implications in terms of 
mixing-enhancement techniques. To explicitly state 
the result, the introduction of disturbances around 
the shear layer in a frequency range anywhere from 
50 to 150 kHz should be able to cause significant 
enhancement of mixing. The frequency range pre- 
scribed here has been obtained after a study of a 
wider range of parameters than discussed herein. 

Figure 13 shows the mean profiles of u versus 
y as a function of axial distance for M c = 0.38. 
The growth of the layer is shown by the dotted 
lines. The distinct change of growth rates with 
distance is noticeable. Since the first part of the 
growth is laminar, a similarity plot on the basis 
of u versus yj\fx is presented in figure 14. The 
laminar similarity seems to be well obeyed for this 
range. A plot containing data up to 100 mm showed 
that no laminar similarity exists in the latter portion 
(not shown here). The turbulent flow similarity is 
on a basis of u versus {y - Vo) / {x — xq). In this 
coordinate, yo is the center of the mixing layer and 
xo is the virtual origin. The value of xq is not known 
a priori. Figure 13 shows that the virtual origin xq 
lies at approximately 0.015 m. Collapsing the data in 
terms of this similarity variable is shown in figure 15. 


Determination of convective speed is based on the 
speed of the structures. The movement of the vorti- 
cal structures also causes a movement of the pressure 
versus x profiles. The plots of p versus x at a spe- 
cific y station for M c = 0.38 and 0.76 are shown in 
figures 16 and 17. Even at M c = 0.38, the changes 
of pressure are, by no means, insignificant; they vary 
between 0.07 and 0.14 MPa. For M r = 0.76, the pres- 
sure variation is between 0.045 and 0.17 MPa. The 
frec-stream pressure is 0.101325 MPa; therefore, the 
pressure variation over the free-stream value is signif- 
icant. The flow has shocks and regions of expansion. 
One of the interesting features of the flow is that the 
regions of shocks occur between structures and over 
them. The flow pattern seen here is different from 
those observed in reference 6 for air at comparable 
values of M c . The shapes of the pressure x plots 
seem to be the same at various times, certainly at 
the times shown in figure 16. On the other hand, the 
shape of the pressure variation for M r = 0.76, shown 
in figure 17, does not remain the same, and there 
arc serious distortions in the shapes for much shorter 
times than the sweep time. Convective speeds ob- 
tained from the structures based on a procedure de- 
scribed previously arc presented in table 4, and the 
results from the standard formula are also shown. 
There seems to be a difference of 100 m/sec in the 
first tanh case, where M c — 0.38. There is a much 
smaller difference in the boundary- layer case; it may 
be good enough for showing the convective speeds. 
For the high values of A/ r , the difference in the two 
results is even larger (^250 m/sec). In one sense, this 
difference is understandable because of the dilatation 
of structures. In another sense, it reflects the possi- 
ble 1 invalidity of the formula at high values of M c . In 
view of this uncertainty, a calculation was made for 
a composition with 30 percent H <2 by mass in the fuel 
stream and with the convective speed evaluated. In 
this case, the density ratio («5.1) was intended to 
verify the nature of the previous results. However, 
even in this case, the difference between the speed of 
the structure and the result from the formula seems 
different by about 100 m/sec. Perhaps the expected 
accuracy of the formula can be no better than 3 to 
4 percent. In the case of M c = 0.76, the result for 
u c is different by 7 percent, and the basic issue of 
dilatation of structures affecting the result remains. 
Hence, the use of M c . in a manner other than desig- 
nating the condition, needs to be examined carefully. 

Figures 18 and 19 show the growth of the layer 
computed from formulas described previously for u, 
p, and H 2 mass fraction. The thicknesses seem to 
decrease after some growth. The fact that growth 
decreases, is, by itself, not very surprising. Forced 
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shear layers show such trends (ref. 19), but typically 
at less than 10 percent. The dips in some quantities 
in figures 18 and 19, particularly the latter, are more 
than this magnitude. Careful examination showed 
that these dips are related to the way that mean pro- 
files of various quantities behave during transition. 
A small kink in the transitioning mean profile can 
cause the presence of a local gradient smaller than 
in the downstream location. This difference causes a 
slightly larger growth rate. It may likewise turn out 
downstream that, although the growth rate is larger, 
there is a lower estimate of the growth. The combi- 
nation of these events produces curves of the shapes 
set'n in figures 18 and 19. The growth rates that 
can be extracted from these plots are discussed sub- 
sequently, along with results of the boundary- layer 
profile. 

Figures 20 to 26 refer to the statistics of fluctuat- 
ing quantities. Figure 20 shows the variation of the 
mis value of a in the field. The top portion shows 
the variation along x at a particular y station. Al- 
though the information contained in the top portion 
is embedded in the two-dimensional plot, it is ex- 
plicitly brought out for a clearer appreciation. At 
larger intensities of fluctuation, the fluctuations de- 
crease over a part of the flow field before the amplifi- 
cation process takes control. The fluctuation profile 
with y is not symmetric. A front is apparently cre- 
ated («5() mm) and is strengthened to a sharp peak 
on the air side. The intensity of fluctuations seems 
to settle down to values of about 6 to 8 percent, ir- 
respective of the initial fluctuation intensity. This is 
one indication of the fact that the subsequent region 
of the flow is becoming independent of the initial dis- 
tribution of the fluctuation. 

Figure 21 shows the variation of v fluctuations. 
The instantaneous values of the v velocities seem 
very large, as large as 400 m/sec, and they have a 
large deviation on either side of zero. These veloci- 
ties cause the flow to show up as large structures on 
a vorticitv plot. The velocity- vector plots in figures 8 
and 10, where the vector is changing directions, are 
also reflective of these cross-stream velocities. The 
changes in velocity directions in a supersonic flow re- 
sult in shock waves, and these cause large pressure 
variations in the mixing layer. The mean velocity, 
however, is not large by comparison (±15 m/sec). 
The peak rms value of a v fluctuation is 8 percent, 
which is comparable to the peak rms value of u fluc- 
tuations. This feature is also reminiscent of fully 
turbulent flows. The fluctuation peaks toward the 
air side (lower). The mean and rms density are pre- 
sented on an x-y plot. The initial change of density 
from 0.075 to 0.175 kg/nr* is evident from the mean 


density plot. The growth of the layer from zero thick- 
ness is clearly visible. The rms value of the p profile 
acquires a top-hat profile with a maximum of 34 per- 
cent. The intensity aspect of the density profile is 
that the top-hat profile is acquired around the tran- 
sition region; broadening of the profile is the only 
thing that happens farther downstream. 

The mean and rms temperature variations are 
shown in figure 23. The temperature profile in the 
present calculation is, like a passive scalar and unlike 
H 2 mass fraction, affected by gas dynamics. Then? 
are changes in the mean temperature to the extent 
of 25 to 30 K, a decrease from the set initial value. 
The extent of the change? of temperature' must be 
understood in terms of energy exchange between 
random energy and orelereel energy. If the inviscid 
global gas dynamics equation is used along stream 
lines as follows: 

l / 2 

c p T -| — — = Constant (14) 

we can differentiate; then, AT ~ —VAV/cp, where 
V is the magnitude of the vector velocity. For typical 
values (V « 2250 m/sec, AV ~ 20 111 /sec, and Cp « 
1800 J/kg-K) AT ~ 25 K, a. value consistent with 
the change shown in the upper plot of figure 23. The 
fluctuation of the rms value of T seems to stabilize 
around 3 percent. This fluctuation is induced by the 
fluid mechanical effects; there are no other energy- 
exchange mechanisms. Figure 24 shows the plot of 
mean and rms values of H 2 mass-fraction profiles; 
rms mass fraction is not normalized. Therefore, the 
peak rms normalized bv the mass fraction of the 
input fuel stream is about, 0.2 percent. While the plot 
of the mean profile looks similar to that of density 
(in a reversed manner), the variation of the rms 
H 2 mass fraction seems quite different from that of 
density. The curve 1 looks like a Gaussian distribution. 
The relationship between mean-density and mass- 
fraction profiles is expected because density changes 
are caused only by molecular weight changes induced 
by H 2 diffusion into the airstream. 

Figure 25 shows the plot of vorticity and Reynolds 
shear-stress fields. The Reynolds shear-stress pro- 
file appears to be growing in size and magnitude at 
100 mm (640 downstream). This trend shows that 
the transition process is incomplete. An attempt is 
made to examine the validity of the Baldwin-Lornax 
model for turbulenc e. In this model, the shear stress 
is represented by pu f v f = Constant 
where Sq is the vorticity thickness, Q m is the mean 
vorticity, arid primes denote fluctuations about the 
mean. The constant is calculated from the actual 
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values of the left- and right-hand terms and plotted 
in figure 26. These calculations were made for both 
101 and 151 grid points to establish the validity of 
the results. As shown in the plot, the "constant is 
constant over a region of the mixing layer but de- 
parts from it significantly, particularly in the peak 
shear-stress area. The fact that the flow has not de- 
veloped into a fully turbulent regime may be partly 
responsible for this departure. If a magnitude must 
be used for the constant, it would be about 0.015. 
Furthermore, what is still not clear is what influence 
the assumption of constancy has on the prediction of 
growth rates of turbulent shear layers. These issues 
have to wait for calculations of modeled equations. 

Boundary-Layer Profile 

Figure 27 shows the plot of the vorticity contours 
for a range of initial disturbances from 0 to 3.68 per- 
cent. The contour plots of vorticity appear thinner 
than for the tanh case, because, unlike the tanh case, 
vorticity changes signs across the shear layer. Note 
the difference in the scales between the boundary- 
layer and tanh cases. 

The plots are very informative. With no im- 
posed fluctuations, there seems to be an onset of 
instability at about 80 mm from the splitter plate, 
and a few vortical structures have begun to develop. 
The spacing of the vortical structures seems much 
smaller than for the tanh case. The average spac- 
ing for the boundary-layer case is about 6 mm and 
23 mm for the tanh case. The spacing is normal- 
ized with respect to initial- layer thickness (2 mm) 
and is 3.0 for the boundary-layer case and ll.o 
for the tanh case. The observed spacing ratio is 
about 2 (ref. 13). This is the first, indication that the 
results of the boundary-layer profile show more real- 
istic results than the tanh profile. The point where 
the first vortical structure gets developed, which may 
be taken as indicative of the beginning of transition, 
moves substantially to lower x values with the impo- 
sition of small, but realistic, disturbance levels. The 
transition distances are 35, 27, 25, and 20 mm at 
progressively increasing disturbance levels. 

The single large jump in the transition distance 
occurs from 80 mm to 35 mm for a change in dis- 
turbance level from 0 to 0.115 percent. Subse- 
quent increases in disturbance level seem to have 
progressively less influence. The transition dis- 
tance, normalized by momentum thickness, is 345 
for a 0.115-percent disturbance level and 265 for a 
0.23-percent disturbance level. The experimentally 
observed transition distance is 165 to 3000 (table 1). 
Transition distances up to 7000 are quoted for various 


cases in reference 13. The transition Reynolds num- 
ber can be computed based on Rp — PwgUc x T / fays- 
When p aVR = 0.125 kg/m 3 , u ( . - 2150 m/sec, and 
XT = 0.04 in, R t is 0.18 x 10 6 . The value of u c used 
here is from the standard equation and is not the 
value obtained from the calculations. The difference 
between the two values is no more than o percent, 
and this difference is quite small for the estimation 
of R. The transition Rp given in reference 9 is 
0.2 x 10 6 ± 0.03 x 10 6 . The results from reference 18 
were quoted previously to be between 0.36 x 10 and 
0.53 x 10 () . These values are based on free-stream 
velocity. If they are based on the convective speed 
of the structures, which, in this case, is half the 
free-stream speed, the transition Reynolds number 
is 0.18 x 10 6 to 0.26 x 10 6 . Again, the comparison 
seems good between the experimental values and the 
predictions. Thus, the present results with boundary 
layer as the initial profile seem to predict realistic 
transitional features at realistic disturbance levels. 

The transition behavior is more vividly depicted 
on time-spectra plots shown in figure 28 for three 
x and y stations. The spectra are shown for both 
0- and 0.115-percent disturbance-level cases. The 
upper plot shows the time spectra at x = 0. The 
disturbance seen in this figure at x — 0 is consistent 
with the description of initial-disturbance structure 
discussed previously. The peak disturbance is about 
1 to 2 m/sec in the central zone. At x — 50 mm down- 
stream, the amplitudes are 10 to 20 m/sec and are 
50 to 90 m/sec at x = 100 mm downstream. The 
mean velocity has also increased from near zero 
to as high as 1500 m/s ec during this transition. 
The nature of the plot at the downstream stations 
shows that most other frequencies, both subharmonic 
and the higher harmonics, are excited perhaps more 
strongly than indicated with the tanh profiles. The 
plots, with respect to the case involving no distur- 
bances, show that the flow has picked up noise that 
was generated numerically. Interestingly, the distur- 
bance is selectively amplified, and peak growth oc- 
curs at about 150 kHz, which implies that the shear 
layer is sensitive to disturbances at this frequency. 
This result can be expressed in terms of a nondimen- 
sional quantity, St = Of /u c , where / is the frequency. 
For 0 = 0.102 mm and u c = 2150 rn/sec, St = 0.007. 
Perhaps St ~ 0.007 as a criterion for excitation is 
valid over a much larger range of parameters. In 
reference 24, a quantity called “preferred mode” is 
indicated in jets and mixing layers, and the Strouhal 
number is between 0.002 and 0.012. It is stated in ref- 
erence 24 that "mixing layers can be manipulated ef- 
fectively with very low forcing levels -0.01 to 0.1 per- 
cent of the average velocity - provided the excitation 
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is applied at the proper frequency.’ 1 This comment is 
entirely consistent with the observations of the 
present computations. 

Figure 29 shows the plots of pressure versus x for 
three successive times. The extent of the pressure 
fluctuations is the same as in the tanh case, but, 
consistent with the vorticity plot, the profile has 
many more waves within the same range of x. The 
convective speed, computed from the movement of 
the pressure waves with time, is presented in table 4. 
As was also indicated previously, this value compares 
well with that from the formula. 

Figure 30 shows the plot of the mean u pro- 
files and the minimum u velocity versus z. A 
plot of minimum u with a nondimensional x de- 
fined by x mn = x/(0Rq) is presented in reference 11 . 

1 he Reynolds number uses the average of the unit 
Reynolds number on both sides. The plot shows that 
at about x mn — 1, the asymptotic value is reached. 
Substituting the values of various quantities and ob- 
taining x, a value of x that corresponds to x mn = 1 , 
yields x — 45.5 mm. The upper plot of figure 30 
shows that the asymptotic value of u seems to have 
been reached at x = 40 to 50 mm. 

Figure 31 shows the plot of the rms value of u 
over the field. The rms value of u increases from 
a relatively low value of 0.1 percent, builds up to 
around 10 percent, and then seems to saturate at 
this level throughout the field after the transition 
process starts. A careful examination shows that 
the rms value of u versus y seems to have attained 
similaiity around 90 mm (900$). The rms value 
of v seems to peak and decay downstream of 600$ 
(fig. 32). The mean v is near zero throughout the 
field. The instantaneous v goes up to ±100 m/sec. 
The mean and rms values of density are shown in 
figure 33. Unlike the tanh case, the rms profiles do 
not seem to have the top-hat shape, and the variation 
with y is asymmetric. The peak rms goes up to 
20 percent, and the mean-density profile appears 
to be approaching similarity. The temperature and 
H ‘2 mass- fraction profiles are presented in figures 34 
and 35. The mean profiles are similar to those 
obtained for the tanh ease. The rms profile shapes 
are slightly different, but the features — namely, the 
peaking of rms intensity and subsequent decay in the 
axial direction and magnitudes arc about the same 
as for the tanh case. 

The growth data for the boundary layer are pre- 
sented in figure 36. As explained previously, the 
growth was obtained by using the u profile. The 
results for most fluctuations show a region of small 
decay before the growth is established; this decay 


may be due to a coupling between the inflow dis- 
turbance structure and the flow instability. For large 
fluctuations, the growth rate is approximately linear. 
In the other cases, a growth rate can be derived in 
the range where there is no decay or constancy. Sur- 
prisingly, the slope obtained for these cases is nearly 
the same. These growth-rate data are presented in 
table 5. The results obtained for other tanh cases are 
also shown in this table. The growth-rate trends for 
various fluctuation intensities are not the same. In 
all cases, the growth rates were estimated from the 
slopes in the latter part of the x range. The linear 
growth rate for the lower intensity fluctuations may 
not have been achieved. The growth rates based on 
H 2 mass fraction or density are the same for all the 
fluctuation intensities. The data for M c = 0.76 were 
derived from similar considerations. Again, there are 
variations. An estimate of the mean value of the 
growth rate was obtained and is included in the ta- 
ble. Further, the growth rate is divided by that for 
the incompressible case (eq. (1)) to obtain the com- 
pressibility effects as in the results of reference 13. 
These data are presented in figure 37 with experi- 
mental results. While there is significant variation of 
the numerically obtained result, the reduction of the 
growth rate with M r is unmistakable. It may well be 
that if the simulations are carried out for much larger 
values of z m /$, the growth-rate data would lead to 
a more consistent picture. 

Figure 38 shows the terms of the conservation 
equations of u, v , and energy. The only dominant 
terms arc the unsteady and convective terms. The 
pressure-gradient terms are about 25 percent of the 
convective terms. The viscous terms seem insignifi- 
cant on the scale shown, but the flow is not necessar- 
ily inviscid. The flow is affected by fluctuations not 
evident in the plot of instantaneous quantities. These 
fluctuations lead to a stress like the viscous stress 
(the Reynolds stress), which is known to be impor- 
tant. To determine their importance, the terms must 
be averaged and the relative orders of magnitude 
must be examined. This determination shows that 
the Reynolds stress terms constitute up to 30 per- 
cent of the convective terms. Figure 39 confirms the 
validity of these results in terms of the variation of 
pressure; variation of pressure affects Reynolds num- 
ber. Changing the Reynolds number by a factor of 4 
causes little change to the vortical structures and the 
mean growth rates. Data for incompressible shear 
layers for local Reynolds number with fixed x (in 
the case considered 0.45 m downstream), in which 
a change in Reynolds number from 4.3 x 10 4 to 
6.7 x 10 4 shows a visible change in growth rate, are 
presented in reference 10. Since the full experimental 
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details of this work are unpublished, the resolution 
of the issues involved in the experiments must await 
discussion at a future time. 

Conclusions 

This study is an exploration of the direct numer- 
ical simulation of high-speed mixing layers. After 
establishing the validity of the approach to resolving 
the mixing layer and the methods of exploring the 
structure of the mixing layer, the issues concerning 
the effects of free-stream disturbance intensity, the 
nature of the profiles, and the role of Reynolds num- 
ber are studied. The following conclusions are drawn 
from the present study: 

1. The calculations show that the profile of the 
hyperbolic tangent calls for unnaturally large distur- 
bance levels («3 to 4 percent) to initiate transition 
in the flow. Boundary-layer profiles, on the other 
hand, show transition Reynolds numbers of about 
0.18 x 10 6 ; these numbers are consistent with the ex- 
perimental results at disturbance levels of the order 
of 0.1 percent. The transition Reynolds number does 
not strongly depend on disturbance intensity over a 
range up to 0.2 percent. 

2. The asymptotic turbulence intensities attained 
in the shear layer for the boundary-layer case (10 per- 
cent) are larger than for the hyperbolic-tangent 
(tanh) profile (6 percent), and the same trend is true 
for Reynolds number stress. 

3. Convective Mach numbers derived from the 
pressure versus streamwise coordinate plots show sig- 
nificant differences from those computed with stan- 
dard formulas to a large extent at high convective 
Mach numbers. These differences are related to the 
significant changes in the shape of the vortical struc- 
tures in the flow field with time. This feature makes 
the idea of a convective Mach number less represen- 
tative in distinguishing flows of different compress- 
ibility at high convective Mach numbers. 

4. The estimates of the growth rates from the 
latter region in growth versus streamwise coordi- 


nate plots are straightforward in most cases, and the 
growth rates that are normalized with incompress- 
ible growth rates decrease with the convective Mach 
number. This decrease is consistent with known ex- 
perimental results. The spread in the results is prob- 
ably the result of the limited axial extent of the field 
explored. 

5. Reynolds number has negligible effects on the 
growth and structure of the mixing layer for the range 
of Reynolds number considered. This result is not 
surprising in view of the dominance of inviscid effects 
over viscous effects. 

6. The utility of the direct simulation approach in 
examining turbulence models is briefly demonstrated 
on the Baldwin-Lomax model. The constant in the 
relationship between the Reynolds shear stress and 
the mean vorticity varies significantly over the region. 
For a subset of the region, the constant remains 
about the same from a value of 0.015. These kinds 
of studies can be performed with the data base that 
is archived for this purpose. 

The relative utility of tanh and boundary-layer 
profiles in studies of mixing layers needs to be con- 
sidered. The tanh profile has the obvious advantage 
of leading to analytical solutions to stability prob- 
lems that are vital for making generalizations. The 
boundary- layer profile, which is closer to reality in 
experiments, has the disadvantage of not being close 
to the solution of the steady problem; therefore, it 
is not amenable to classical stability treatment. If 
one intends to capture the transition processes and 
compare them with experiments, there is no alterna- 
tive to the boundary-layer profile. Thus, the tanh 
and boundary-layer profiles both have their regions 
of importance in the study of mixing layers. 


NASA Langley Research Center 
Hampton, VA 23665-5225 
April 1, 1992 
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Appendix 

Numerical Algorithm 

For supersonic Hows, no other algorithm has been 
more widely used than the method discussed in ref- 
erence 25. Over the last 20 years, this method has 
been repeatedly compared with other codes and with 
experimental data. Recent efforts have been made to 
improve the accuracy of the method by increasing the 
spatial truncation error to fourth order and thus re- 
ducing some of the dispersive errors of the method. 
A fourth-order variant of the original MacCormack 
algorithm (a five-point linear stencil) has been ad- 
dressed in references 22 and 26. These references 
show that increased resolution is obtained for a broad 
class of problems with the fourth-order method. 

In the last decade, a wealth of upwind algorithms 
have been developed. These methods incorporate nu- 
merical stencils which align with the direction that 
information is propagating in the flow. This prop- 
agation is usually accomplished by decomposing the 
governing equations into a form in which each charac- 
teristic wave can be differenced in an upwind manner. 
Numerous implementations of these methods have 
shown a great deal of success in resolving supersonic 
problems of practical interest. In spite of the success 
of upwind algorithms for the Euler equations, difficul- 
ties exist for their use in supersonic chemically react- 
ing flow problems. Specifically, the chemistry system 
is extremely difficult to formulate along characteristic 
directions, and many of the upwind methods are not 
easily generalizable with chemically reacting flows. 

The basic characteristics of central- and upwind- 
difference algorithms, as well as the beneficial at- 
tributes of higher order algorithms, can be observed 
in a two-dimensional model problem (referred to as 
the color problem) designed to demonstrate numeri- 
cal diffusion and dispersion. (See ref. 27.) Different 
numerical algorithms are compared with an exact so- 
lution and with each other to determine the charac- 
teristic of each. The test problem used in this study 
is the two-dimensional advection of a scalar field de- 
fined by the equation 


with 


d 4 > 
~dt 


djuv) d(v*) 

Ox + dy 


^(x,i/,()) = *o{x,y) 


(15) 


(16) 


and with suitable boundary conditions. The velocity 
vector was chosen to be a solid-body rotation defined 


by 


U = -Sly 1 
V = Six ) 


(17) 


with Q = 27 r, so that one complete revolution of the 
flow occurred in one unit of time. The domain that 
was used corresponded to the domain — 1 < x, y < 1. 
Dirichlet or first-order extrapolations were used for 
the boundary conditions for flows entering or exiting 
the domain. The boundaries were located suitably 
far from the regions of interest, so that the lower or- 
der treatment did not affect the quality of the results. 
The grid used 100 by 100 uniformly spaced points on 
which to solve the governing finite-difference equa- 
tions. All calculations were integrated in time un- 
til the initial distribution had rotated exactly once. 
The distribution was then compared with the initial 
distribution by using two error norms: the L ^ or 
maximum error norm, and the L 2 norm. 

Figure 40 shows the initial distribution on the 
grid used in this study. The initial value of 4> was 
chosen to be 1 in the interior of the “L-shaped” body 
and 0 elsewhere. This distribution is discontinuous 
and is a serious test of the dispersive properties of an 
algorithm. Reference 28 shows that this discontinuity 
seriously degrades the accuracy of the solutions and 
that the error (A x or L 2 ) decays no better than first 
order with increasing spatial resolution, regardless of 
the spatial accuracy. With higher order algorithms, 
the error is confined to a narrower region about the 
discontinuity and thus contributes to a slightly better 
L ‘2 norm. This distribution was chosen to more 
closely approximate the discontinuous behavior of 
the Euler equations. 

Figures 41 and 42 show the results obtained after 
one complete revolution of the initial distribution 
with the standard MacCormack and fourth-order 
Gottlieb- Turkel algorithms. Contours ranging from 
0 to 1 arc plotted in steps of 0.05 (21 contour levels). 
The sharp distinction between regions initially at a 
value of — 1 and the surrounding regions at = 0 
is no longer apparent. The dispersive nature of the 
algorithm has distorted the solution dramatically. 
The regions inside the innermost contour level are 
at values greater than 1 and are numerical in nature. 
The convection of the initial discontinuity is much 
better accomplished with the fourth-order method. 
The outer boundary of the initial distribution is 
nearly replicated with the fourth-order schemes. The 
contour levels in the interior regions of the solution, 
where 4* < 1, indicate that the numerical dispersion 
is present. Thus, although the fourth-order schemes 
track the discontinuity with much greater accuracy, 
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they still exhibit significant over and undershoots 
near the discontinuity. 

Figure 43 shows the same test problem with the 
third upwind-biased algorithm. In all these cal- 
culations, the stencil was oriented in the upwind 
direction, and was therefore varied in each of the 
four quadrants of the test problem. The third-order 
upwind-biased algorithm produces a clean solution 
after one complete rotation. The location of the dis- 
continuity is somewhat diffused at the corners rela- 
tive to the fourth-order central-difference algorithms, 
but the two-dimensional dispersion present in the 
central-difference algorithms is eliminated. 

When the quantitative features of the upwind al- 
gorithms are compared with the central-difference al- 
gorithms, it is apparent that the fourth-order scheme 
tracks the discontinuities more accurately than does 
the upwind algorithm, but it has far more two- 
dimensional overshoots and undershoots than the up- 
wind algorithm. Except for the diffused corners of 
the discontinuity, the third-order upwind algorithm 
provides a very accurate solution. The algorithms 
of at least third-order accuracy provide significant 
improvement over the second-order central-difference 
algorithms. 

To quantify some of these observations, the error 
after one complete revolution of the initial distribu- 
tion was calculated. Both the norm (shown by 
the maximum or minimum of the solution over the 
entire domain) and the L 2 norm arc presented. Ta- 
ble A1 shows all the algorithms, the maximum and 
minimum values of 4 / , the L 2 norm, and the CFL 
at which the calculations were performed. The L max 
and L m - m norms for the upwind case are reduced by 
the higher order algorithms, especially for the upwind 
algorithms. Also, the L 2 norm for the third-order 
upwind and the fourth-order schemes is nearly the 
same. These observations arc characteristic of the 
comparison between upwind and centered schemes. 
With higher order schemes, the efficiency and robust- 
ness achieved in solving realistic problems decrease. 
The third-order upwind-difference and fourth-order 
central-difference schemes represent truncation lev- 
els that compromise the extremes of accuracy versus 
efficiency. 

The SPARK2D code used in the present study 
incorporates a combination of upwind- and central- 
difference schemes. The numerical simulation of the 
supersonic mixing layer is ideally suited for the use 
of upwind-biased algorithms in the streamwise direc- 
tion. Because the flow is supersonic, all the eigenval- 
ues are oriented in the same direction. The equations 
can be differenced in an upwind manner without a 


characteristic decomposition. A third-order upwind- 
biased algorithm has been chosen as the numerical 
algorithm for the streamwise direction. It shows low 
phase errors and acceptable levels of numerical dissi- 
pation, with extremely good characteristics for cap- 
turing discontinuity. 


Table Al. Error Norms From Different Schemes 


Scheme 

Order 

I-'inax 

^min 

1-2 

CFL 

MacCormack 

Af~. As* .Ay 2 

1.429 

-0.378 

7.7(> x nr 4 

1 

(ref. 25) 
Gottlieb- Turkel 

At 2 . Ax', Ay' 

1.337 

-.239 

5.19 

2/3 

(ref. 22) 
Third upwind 

At', Ar', Ay' 

1.148 

-.108 

5.29 

1 


Extremely strong gradients exist in the cross- 
stream direction. A higher order numerical algorithm 
is ideally suited for this direction. Since the flow in 
this direction is always subsonic, central-difference 
methods are ideally suited. The fourth-order method 
of reference 22 is used in the cross-stream (and nor- 
mal for the three-dimensional case) direction. The 
numerical algorithm in two spatial dimensions can be 
described by the model equation Uf, + F x T Gy = 0 
as 


lJ ij = u lj ~ Jxz ( 2f '+h; + * F '-J ~ 6F '- l '.y + 


l 'J 6 Ax 

At 


6 Ay 




At 


t/r . = _[^ + ^__(2F i+lj+ 3F, 


-6Fj_i j + Fi- 2j) 




(18) 


The fourth-order Gottlicb-Turkel scheme is stable for 
CFL < | and thus limits the maximum allowable 
CFL. The Gotti ieb-Turkel scheme is not applicable at 
the grid point next to the boundary and must rely on 
a standard MacCormack difference formula at that 
point. Similarly, the upwind scheme is not defined 
at the first and second axial grid points and relies on 
a first-order upwind scheme at that point. Explicit 
numerical damping, which is often added to codes 
for engineering calculations, was not used here for 
reasons of accuracy. The only numerical dissipation 
present was that generated by the algorithm itself. 
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Code Validation 

The first step in validating the code was to deter- 
mine if the numerical methods used in the SPARK2D 
code are capable of resolving the linear growth of 
the two-dimensional mixing layer. It is essential 
that all the aspects of the linear regime of the 
mixing layer be well resolved for grids comparable 
to those used in this study before addressing the 
far more difficult nonlinear problem of vortex roll- 
up and chemical reaction. Linear stability theory 
predicts that the temporally developing compress- 
ible two-dimensional mixing layer (air into air) is 
unstable for a velocity profile that is initially speci- 
fied as a hyperbolic-tangent axial-velocity distribu- 
tion. In the initial stages of the instability (the 
linear regime), the growth of the unstable modes is 
exponential. For these calculations, accurate eigen- 
modes are provided from a spectral linear stability 
code developed in reference 3. From these eigen- 
modes, nondimensional growth rates and characteris- 
tic frequencies were calculated. This nondimensional 
growth rate provided a reliable measure of the accu- 
racy of the finite-difference algorithm being tested. 

The temporal growth of the two-dimensional mix- 
ing layer is simulated by assuming the flow to be pe- 
riodic in the streamwise direction. The period corre- 
sponds to a wavelength of x — 0.6283 mm (the mode 
that grows most rapidly, as determined from the lin- 
ear analysis). The initial velocity distribution is spec- 
ified as U{x,y, 0) = Uoc tanh (|), T(x,y,0) = 0.0, 
T{x,y,0) ~ 293 K, and P(x,yA)) = 1 atm. Source 
terms were added to the Navier-Stokes equations, so 
that the momentum and energy equations would pre- 
serve the free stream. The resulting Reynolds num- 
ber of the flow is 0.187 x 10'*, based on the layer thick- 
ness, where U 3 c = 100 m/sec. The Mach number is 
0.30. The width of the half-layer D is 2.5 x 10~ 2 mm, 
and the half- width of the domain is 100 times the 
layer thickness. The grid in the streamwise direction 
is uniform, while the grid in the cross-stream direc- 
tion is highly stretched. The transformation of the 
grid ensures that about half of all the y grid points 
are located within the initial mixing-layer width. 

For these studies, the unstable modes in the layer 
were allowed to grow from the numerical instabilities 


produced by the machine round-off errors. Expo- 
nential growth of these modes was seen in all cases 
after an initial transient period. The “linear” regime 
was characterized by the growth period, during which 
the product of perturbation quantities was still neg- 
ligible. For these tests, that period was arbitrarily 
between 2.5 x 10 0 and 2.75 x 10 - ^ sec and corre- 
sponded to disturbance amplitudes three to five or- 
ders of magnitude larger than machine round-off er- 
rors. A grid convergence history of the numerical 
method was used to determine its formal accuracy 
and to ascertain the grids necessary to resolve this 
fundamental phenomenon . 

The most unstable mode in this problem grows 
at an exponential rate, with an exponent determined 
from linear stability theory to be 0.140000 in non- 
dimensional units. This rate was used as the “ex- 
act” growth rate for these conditions. A series of four 
grids were then defined, each with a grid density that 
was a constant multiple of the previous grid. Three 
algorithms MacCormack (ref. 25) at CFL = 1.0, 
Gottlieb-Turkel (ref. 22) at CFL = 0.5, and DCPS 
(ref. 27) at CFL = 0.5 — were then run on identi- 
cal grids, and the nondimensional amplification rate 
was determined from an integration of the energy 
spectrum. (The integral of the fundamental mode 
over the entire domain was monitored in time. The 
change, with respect to time, yields the amplification 
rate.) 

In figure 44, the amplification rates of the three 
methods are plotted against to show quar- 

tic accuracy. The symbols N and AT max are the num- 
ber of grid points and the maximum number of grid 
points used in the study, respectively. Here, the 
fourth-order spatial accuracy of the spatial derivative 
from each- method is indicated by the linear conver- 
gence of the solution to the predetermined amplifica- 
tion rate. The Gottlieb-Turkel scheme and the DCPS 
algorithm both converge with fourth-order accuracy. 
It is apparent from these results that, even for coarse 
grids (10 grid points per wavelength), the fundamen- 
tal features of the linear growth of the mixing layer 
are resolved. There should be no question that the 
grids used in these studies are sufficient to resolve the 
linear regime of the two-dimensional mixing layer. 
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Table 1. Data and Results of Previous Investigations 




Quantity 

Reference 11 
(case 1) 

Reference 11 
(case 2) 

Reference 13 

Reference 16 

Ml /Ml 
Vi/U 2 , m/sec 
61 / 62 , mm 
0l/0*2, 111 m 
P \/ P 2 , Pa 

Ti/n, k 

Me 

2.8/2.15 
655/580 
2. 0/2.0 
0.123/0.152 
2089/2294 
130/174 
0.325 

2.7/1.67 

630/515 

1.6/1. 4 

0.12/0.125 

2895/3673 

134/206 

0.515 

3. 1/1. 7 and others 
Ratio =1.4 
0.44/0.36 
0.042/0.069 
NA 
NA 
0.6 

1.95/1.33 

500/380 
2. 0/2. 3 
NA 

Nearly same 

169/220 

0.31 

Gases 

Air/air 

Air/air 

N 2 /N 2 and other gases 

Air/air 

Comments 

u profile takes more 
than 80 mm to go 
near equilibrium; 
Transition « 1650; 
dominant frequency 
~ 93 kHz 

u profile takes more 
than 80 mm to go 
near equilibrium; 
Transition ^ 3000; 
dominant frequency 
« 130 kHz 

Growth rate by pitot 
surveys; Transition 
« 40 to 60 mm 

Growth rates 
measured 
as 0.011 


Table 2. Inflow Parameters 

[T = 2000 K; p — 0.101325 MPa; 6 (tanh case) — 1.54 mm; 0 (boundary- layer case) = 0.102 mm] 


Composition 

Parameter 

Case 1 

Case 2 

Fuel: 

P , kg/ in 3 

0.075 

0.075 

0.1H 2 + 0.9N-2 

it, m /sec 

2670.0 

2670.0 


M 

2.0 

2.0 


7 

1.3133 

1.3133 


a, m/sec 

1336.0 

1336.0 


6 , 111 m 

1.0 

1.0 


/x, kg /in- sec 

5.5 X 10~ 5 

5.5 x 10 -5 

Oxidant: 

Pi kg/ m :! 

0.175 

0.175 

0.2320-2 + 0.768N 2 

u, m/sec 

1814.0 

4330.0 


M 

2.1 

5.25 


i 

7 

1.296 

1.296 


a, m/sec 

864.3 

864.3 


6 , mm 

1.0 

1.0 


/i T kg/m-scc 

6.2 x 10“ 5 

6.2 x 10 -5 


u c , m/sec 

2150.0 

3670.0 


M c 

0.385 

0.76 
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Table 3. Nondimensional Parameters 


Parameter 

Case 1 

%m/ @ 

64 (100 mm, tanh) 
100 (150 mm, tanh) 
130 (200 mm, tanh) 
1000 (100 mm, b.l.) 
1500 (150 mm, b.l.) 

Temperature ratio (fuel /oxidant) 

1.0 

Pressure ratio (fuel/ oxidant) 

1.0 

Density ratio (fuel/oxidant) 

0.429 

Velocity ratio (fuel/oxidant) 

1.45 

Momentum ratio (fuel/oxidant) 

0.90 

R<5, fuel 

3 640 

R^air 

5120 

Rg = p a u c 0/n, tanh 

5518 

Ra = PnUr0/a, b.l. 

360 


Case 2 


64 (100 mm, tanh) 
100 (150 mm, tanh) 
130 (200 mm, tanh) 
1000 (100 mm, b.l.) 
1500 (150 mm, b.l.) 


1.0 

1.0 

0.429 

0.616 

0.163 


3 640 
12 220 


11 100 
630 


Table 4. Convective Speeds 


Me 

Composition 

u f: , m/sec 

from formula 

u c , m /sec from 
vortical structures 

0.38 tanh 

0.1H 2 ± 0.9N 2 

2150 

2260 ± 50 

0.38, b.l. 

0.1H 2 + 0.9N 2 

2150 

2190 ± 30 

0.76, tanh 

0.1H 2 + 0.9N 2 

3670 

3400 ± 150 

0.38, tanh 

0.3H 2 + 0.7N 2 

2543 

2648 ± 50 


Table 5. Growth Rates 


M r 

Fluctuation, 

percent 

On 

<^ h 2 ’ 

Mean 6 1 



3.68 

0.03 ± 0.003 

0.05 ± 0.002 



0.38, 

1.84 

.05 ± 0.002 

.05 ± 0.002 

0.05 

0.77 to 0.83 

tanh 

.92 

.06 ± 0.020 

.05 ± 0.002 




3458 i 

0.075 ±0.005 

0.07 ± 0.005 



0.76, 

1.84 

.04 ±0.003 

.04 ± 0.003 

0.04 

0.45 to 0.50 

tanh 

.92 

.04 ±0.002 

.04 ± 0.005 




3.68 

0.055 ± 0.002 




0.38, 

.92 

.055 ± 0.005 


0.055 

0.84 to 0.91 

b.l. 

.46 

.055 ± 0.005 





.115 

.055 ± 0.005 








Air 












Amplitude Amplitude Amplitude 


Grid Domain 

201 by 101 100 mm by 30 mm 

301 by 101 100 mm by 30 mm 



Frequency, Hz Frequency, Hz Frequency, Hz 



Frequency, Hz Frequency, Hz Frequency, Hz 



Frequency, Hz Frequency, Hz Frequency, Hz 

Figure 6. Plot of time spectra for u at three x and y stations for control and finer x grid cases. Amplitude is 
the square root of the energy of fluctuations. 
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Amplitude 
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Amplitude 














Velocity vectors 



Figure 10. Vector and contour plots of velocity, vorticity, pressure, and H 2 mass fraction (tanh profile) 
at 3.85 percent fluctuation and M c = 0.76. Grid: 251 by 101; domain: 125 mm by 30 mm. 







t r 


Grid: 201 by 101 
Domain: 100 mm by 30 mm 



Figure 11. Contour plots of vorticity at various i 
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0 percent 
fluctuation 


0.92 percent 
fluctuation 



1 .84 percent 
fluctuation 


T 


T 



t r 



3.68 percent 
fluctuation 


[-disturbance intensities and M c — 0.76. 
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Figure 12. Effect of location of disturbances on vorticity at x = 0 on shear growth layer at 3.85 percent 
fluctuation and M c = 0.38. Grid: 201 by 101; domain: 100 mm by 30 mm. 
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Initial growth 
Asymptotic growth 



Figure 13. Mean profiles of u versus y with axial distance at M , 










.10 



Figure 15. Profiles of figure 13 on a turbulent similarity plot for x = 80 to 100 mm. Domain: 100 mm by 
30 mm; x station (x s t a ) — 16 mm; = 15 mm. 
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-500 



Figure 21. Instantaneous v at three fixed values of y and mean rms v as a function of x and y for 3.68 percent 
fluctuation and M c = 0.38. Grid: 201 by 101; domain: 100 mm by 30 mm; initial profile: tanh. 
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0 percent fluctuation 


t r 1 1 1 1 r 



0.1 1 5 percent fluctuation 



0.46 percent fluctuation 




0.92 percent fluctuation 


t 1 1 1 1 r 





3.68 percent fluctuation 


i j 1 1 1 1 1 r 


Vorticity contour plot for boundary-layer case at various fluctuation intensities and M r = 






Amplitude Amplitude Amplitude 














































Terms 

Figure 38. Magnitude of terms of u, v, and energy equations for boundary-layer case at M c — 0.38. Domain: 
150 mm by 30 mm; 0.115 percent fluctuation. 
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Figure 40. Initial distribution for “color” problem, superimposed on grid used in study. 



Figure 41. Rotation of initial distribution with MacCormack (ref. 25) algorithm. Accuracy of method is 
second-order in time and space. 







0(At 2 ,Ax 4 ,Ay 4 ) 



Figure 42. Rotation of initial distribution with Gottlieb-Turkel (ref. 22) algorithm. Accuracy of method 
second-order in time and fourth-order in space. 



Figure 43. Rotation of initial distribution with third-order, upwind-biased algorithm. Accuracy of method 
second-order in time and third-order in space. 





Figure 44. 
scheme. 
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Comparison of higher order Gottlieb-Turkel and DCPS schemes with second-order MacCormack 
ft* is st ream wise grid density; N y is cross-section grid density. 
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